Resonance induced by repulsive interactions in a model of globally-coupled bistable 

systems 
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We show the existence of a competition-induced resonance effect for a generic globally coupled 
bistable system. In particular, we demonstrate that the response of the macroscopic variable to an 
external signal is optimal for a particular proportion of repulsive links. Furthermore, we show that 
a resonance also occurs for other system parameters, like the coupling strength and the number of 
elements. We relate this resonance to the appearance of a multistable region, and we predict the 
location of the resonance peaks, by a simple spectral analysis of the Laplacian matrix. 
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I. INTRODUCTION 

The amplification of an external forcing acting upon a 
dynamical system under the presence of the right amount 
of disorder has attracted much attention in the last 
decades. The phenomenon of stochastic resonance, pro- 
posed in 1981 to explain the periodicity of ice ages, [l|, 0] 
is a somehow counterintuitive effect arising from the co- 
operation between deterministic dynamics and dynami- 
cal disorder or noise. By this effect, a system's coherent 
response to a weak signal can be optimally amplified by 
an intermediate level of noise. The prototypical example 
is that of a continuous variable whose deterministic dy- 
namics is relaxational in a double-well potential. Noise 
induces jumps between the wells with a rate given by 
Kramers' expression Q. The system becomes optimally 
synchronized with the signal when the signal half-period 
matches Kramers' rate, as reflected by a maximum value 
in a suitably defined response. Applications of stochas- 
tic resonance wer e addressed in many areas, and the 
theory evolved in several directions (see [i| and Q for re- 
views). We would like to stress here the focus on many- 
component, or extended, systems [(| 0] for which it was 
shown that it is possible to tune some of the param- 
eters, like the number of elements [8T4TT1]. or coupling 
strength p^ - |20j . to further enhance the resonance effect. 

Another, more recent, related line of research considers 
the role that other types of disorder, such as quenched 
noise (identified with heterogeneity or disorder), can play 
in producing a resonance effect in many-component sys- 
tems. Tessone et al. [U 22] have shown that in generic 
bistable or excitable systems, an intermediate level of 
diversity in the individual units can enhance the global 
response to a weak signal. In bistable systems, a mean 
field analysis interprets the phenomenon macroscopically 
as a result of the two equilibrium states getting closer to 
each other and the lowering of the height of the potential 
barrier, thereby turning the signal suprathreshold. The 
resonance appears close to an order-disorder transition, 



and it was suggested that any source of disorder would 
lead to the same result. This diversity-induced resonance 
effect has already found extensions in systems concerned 
with cellular signaling [23[ , complex networks , inter- 
cellular Ca 2+ wave propagation induced by cellular diver- 
sity [111, linear oscillators (2(| and others. Focusing on 
the double-well model, Perc et al. studied the com- 
bined effect of dynamic and static disorder, where static 
disorder was either diversity, the presence of competitive 
interactions, or a random field. Namely, they showed 
that the random presence of repulsive bonds decreases 
the level of noise warranting the optimal response. 

It is the purpose of this work to show that competitive 
interactions can actually replace -not merely enhance- 
noise in its constructive effect. In a recent study [28[ 
we focused on a network of Ising-like units, and found 
that the addition of an intermediate fraction of repul- 
sive links can increase the sensitivity to an external forc- 
ing. In this work, we expand the result by looking at 
a generic double-well model, and clarifying some of the 
previous conclusions. We resort to a spectral analysis 
of the Laplacian [29| matrix, to locate the amplification 
region, and unveil the mechanism of resonance. 

The outline of this paper is as follows: in section [TT] 
we will introduce the model; we show that there is an 
amplification and discuss how the amplification mecha- 
nism is related to a break of stability in section IIII1 and 
how we can predict the resonance peaks in section IIVI 
Conclusions are drawn in section \V\ 



II. THE MODEL 

We consider a system of N globally-coupled bistable 
units described by real variables Si(t),i — 1, . . . , N under 
the influence of a periodic forcing. 
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where t is the dimensionless time, C measures the 
coupling strength amongst the different units and 
Asin(27rf/T) is a periodic external signal with amplitude 
A and period T. The interaction matrix Jij reflects the 
presence of attractive and repulsive interactions between 
the units. More specifically, we adopt the following val- 
ues at random: 



— 1, with probability p, 
1, with probability 1 — p. 



(2) 



The single-element case, N = 1, with added noise is 
the prototypical double-well potential system for which 
stochastic resonance was first considered. The case with- 
out repulsive interactions, p = 0, can still be described 
globally by a bistable potential (see next section) and, in 
the presence of noise, has been widely studied as a model 
case for stochastic resonance in extended systems [7J; it 
has also been considered, in the presence of a random 
field, as a prototypical example for the diversity-induced 



resonance effect [21] . For p > 0, the coexistence of at- 



tractive and repulsive interactions is characteristic of a 
wide class of spin-glass- type systems |30j ] . 

We will focus on the macroscopic variable S(t) — 
7? Si s i an d use as a measure of response the spec- 
tral power amplification factor [3lj . defined as the ratio 
of the output to input power at the corresponding driving 
frequency: 



R = AA- 



S{t)) 



(3) 



where (• • • ) is a time average. R is roughly proportional 
to the amplitude of the oscillations of Sit). If R < 1, 
then the amplitude of the response is less than that of 
the signal, and vice versa for R > 1. 



III. SIGNAL AMPLIFICATION 

It is convenient to analyze first the structure of the 

steady-state solutions for the system of equations (J]) in 

the non- forced case, A = 0. The dynamics is relaxational 

dsi dV r— i . 

[32JJ, being 
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the Lyapunov potential. Therefore, the stable steady 
states are the configurations (si, . . . , sjv) which are abso- 
lute minima of V. If there are no repulsive links, p = 0, 
the Lyapunov potential has just two equivalent minima 
at Si = +1 or Si = —1, V« = 1,...,N and, hence, the 
macroscopic variable will reach the stable asymptotic val- 
ues 5 = +lor5=— 1, a typical situation of bistability. 
Whether S = +1 or S = — lis reached depends exclu- 
sively on the initial conditions. As p increases, the abso- 
lute minima depart from S — ±1 and, furthermore, new 
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FIG. 1: We plot in a gray scale the stationary probability 
distribution -P(S'), in the absence of external signal A = 0, 
coming from numerical simulations of Eqs. fTJ. For better 
viewing, the distribution has been rescaled by its maximum 
value at each p. The data show that at p < p c ~ 0.44 the sys- 
tem presents two equivalent absolute maxima for P(s), while 
there is only one absolute maximum for p > p c . We note, 
however, that there are many relative maxima for all values 
of p, specially around the region p ~ p c . Other parameter 
values are: N = 200, C = 8. The probability has been com- 
puted after averaging over 1000 realizations of the couplings 
Jij and initial conditions drawn from an uniform distribu- 
tion in the interval ( — 1, 1). For the numerical integration we 
used a fourth-order Runge-Kutta method with a time step 
At = 0.1. 



metastable minima of V appear. The dynamical equa- 
tions ([l} may or may not get stuck in one of these min- 
ima, depending on initial conditions and the particular 
realization of the coupling constants Jy. We have used 
throughout the paper random initial conditions drawn 
from a uniform distribution in the (—1, 1) interval, al- 
though we have observed the same type of phenomenol- 
ogy when using other random , but still symmetric, dis- 
tributions such as truncated Gaussian or the Johnson 
family of distributions. 

/From our simulations we compute numerically the 
probability distribution P{S) of the final values of S 
reached during the dynamical evolution for different re- 
alizations of the coupling constants Jy and initial condi- 
tions. This is plotted in Fig.[TJ We can observe a second- 
order phase transition as the average value (IS'(t)l) van- 
ishes for p > p c w 0.44. One can interpret these results in 
terms of an effective potential V c s(S) = — \n.P(S) which 
has two equivalent absolute minima at S = ±So{p), 
where 1 > So(p) > for < p < p c , and one absolute 
minimum at S — for p > p c . The effective potential 
VcS presents many relative minima for all values of p > 0, 
especially in the critical region p w p c , a typical situation 
for the spin-glass models [30|. 
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FIG. 2: Spectral amplification factor 7? versus probability of 
repulsive links p. Inset: the influence of the amplitude A of 
the external forcing on the response R. For suprathreshold 
amplitudes, A > 0.4, R decreases with A due to the denomi- 
nator A 2 in the definition of the spectral amplification factor 
R. T = 300, N and C as in Fig. □ 



We now turn on the forcing A > and study the sys- 
tem response, as measured by the spectral power ampli- 
fication factor R defined above, Eq. ([3]). Consider first 
the case p = 0. For a small, sub-threshold, amplitude 
A the macroscopic variable S(t) will just execute small 
oscillations of amplitude proportional to A around the 
stable values S = +1 or S = — 1. As A increases beyond 
the threshold value A w 0.4 the amplitude of the forcing 
is large enough to induce large jumps of the macroscopic 
variable from S ~ — 1 to 5 ~ +1 and vice versa. This 
change of behavior at A a appears as a sudden increase 
in the value of R, as shown in the inset of Fig. [2] As 
the same inset shows, similar behavior is observed for 
< p < p c : the response shows a sudden increase for a 
particular value of the amplitude A and then decreases 
monotonically. For p > p c , the response is very small and 
almost independent on the value of A. 

More interesting, and the main result of this paper, 
is the dependence of R on the probability p of repulsive 
links, main plot in Fig. [2] We note that there is an opti- 
mal probability of repulsive links that is able to amplify 
signals whose amplitude would be sub-threshold in the 
case p = 0, i.e. A < A a . For suprathreshold signals, 
A > j4 g , the presence of repulsive links does no longer 
lead to enhanced amplification. As shown in the figure, 
the optimal value for amplification is close to the criti- 
cal value p c signaling the transition from bistability to 
monostability in the non-forced case. The optimal am- 
plification as a function of p can clearly be observed in 
Fig.® which shows representative trajectories for p = 
(small oscillations around the value S = +1), p = p c 
(large oscillations between S » +1 and S » —I) and 
p = 1 (small oscillations around 5 = 0). 

The existence of an optimal value of the fraction of 
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FIG. 3: Representative trajectories of the macroscopic vari- 
able S(t). Note the large amplitude of the oscillations in the 
intermediate case p = 0.44. The "signal" is the periodic func- 
tion Asin(2Trt/T). Values of N, C and T as in Fig.d A = 0.2. 
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FIG. 4: The influence of the signal period T on the response 
R. In the inset, we see the response R reaches a constant 
value for slow enough signals. Values of N, C and A as in 
Fig. El 



repulsive p for which signal amplification is maximum 
is somehow reminiscent of the stochastic resonance phe- 
nomenon. There are some important differences, how- 
ever. While in stochastic resonance, the response R 
shows a maximum as a function of period, resulting from 
the matching between Kramers' rate and the forcing half- 
period, in our case the same optimal disorder p amplifies 
responses to signals of every period, as shown in Fig. U 
When the signal is slow enough, the system has time 
to respond to the fuller extent, going to the absolute 
extrema of the potential, and the amplification factor 
reaches a constant value, see inset of Fig. @] 

It is possible to reinterpret these results in terms of 
the effective potential V e g(S) introduced above. The pe- 
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riodic forcing can be seen, approximately, as a periodic 
modulation to the potential V c $ — S ■ A sin(27rf/T). As 
discussed above, the effect of the repulsive links is such 
that V e s(S) changes from bistable at p = to having 
many metastable minima at p « p c and a single absolute 
minimum for p > p c . Hence, the deep potential barrier 
separating the S = ±1 solutions for p = lowers under 
the effect of the repulsive links. As a consequence, the 
modulation induced by the periodic forcing is now large 
enough, and the global variable is then able to oscillate 
from the minimum +Sg(p) to —Sq(p) and vice versa. As 
p approaches p c a more complicated scenario appears. 
In this region the effective potential presents already a 
rich structure with many metastable minima in the non- 
forced case. Those minima can be mod ified or even 
disappear by the effect of the periodic modulation. It 
is particularly illustrative to compare the responses to a 
suprathreshold signal of amplitude A = 0.4 in the case 
p = 0, and to a signal of amplitude A — 0.2 (which would 
be subthreshold in the case p = 0) at the optimal fraction 
of repulsive links p — p c . In both cases, the amplitude of 
the oscillations is approximately the same, as the system 
makes large excursions from S ~ — I to 5* « +1 and vice 
versa. However, the shape of the oscillations is rather 
different, as shown in the upper panel of Fig [5j In the 
p = case, the transition from one minimum to the other 
is rather fast (vertical portion of the dashed line), while 
in the case p — p c , the transition is slower as the system 
seems to be spending more time in intermediate states. 

To determine what those differences reveal about the 
underlying effective potential, we have used a method [33[ 
that allows us to detect the number of states a system 
visits from an analysis of its time series. A typical ex- 
ample is shown in the lower panel of Fig. [SJ We only 
detect two states in the global variable S when p = 0, 
corresponding, as expected, to the modulated bistable 
potential. By contrast, the slight irregularities in the tra- 
jectory for p — 0.44, hardly visible by eye, correspond to 
several very shallow potential wells. The system evolves 
through many states at the optimal probability of repul- 
sive links, as shown in the lower panel of Fig [5J This 
image explains why signals of every amplitude and pe- 
riod can be amplified for p w p c . In this case, the system 
can access the many intermediate states, covering a dis- 
tance proportional to T and A, in case of very fast or 
very weak signals. 

The previous results show that the disorder induced 
by an intermediate level of repulsive links is an essential 
ingredient to get an optimal response to the external forc- 
ing. This can be explained as, in the absence of forcing, 
the metastable states correspond to a wide distribution 
for the values of s^s. When the forcing is turned on, some 
units will be responsive to the signal, and then they will 
pull others which are positively coupled to them. This 
basic mechanism is further highlighted by the observation 
of a resonance behavior with both the coupling constant 
C and the number of units N. 

The resonance with C and some representative trajec- 
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FIG. 5: We amplify some representative trajectories (upper 
panel), and count the number of states through which the 
system moves in each trajectory (lower panel). Values of T, 
N and C as in Fig. H 
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FIG. 6: Coupling-induced resonance. Main plot: the response 
R shows a maximum as a function of coupling constant C. As 
shown in the insert, the same maximum appears as a function 
of the localization measure M, see section IIVI Values of N, 
T and A as in figure [3] and K — 0.2 in the insert. 



tories are displayed in Figs. |5] and [71 respectively. In the 
weak coupling limit, the units behave basically as inde- 
pendent from each other and, as the signal amplitude A 
is subthreshold for a single variable, the overall response 
is small. In the large coupling limit, the interaction term 
is too big to allow an unit that could first follow the signal 
to depart from the influence of its neighbors. 

The resonance with the number of units N and some 
representative trajectories is presented in Figs. [5] and M 
Since fluctuations in the number of repulsive links de- 
crease with N, a larger system requires a greater fraction 
of repulsive links to achieve the same level of disorder 
than a smaller system. As a consequence, the response 
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FIG. 7: Coupling-induced resonance, as revealed by the reso- 
nant trajectory at optimal G = 11. Values of N, A and T as 
in Fig. H 
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FIG. 8: System-size induced resonance. Main plot: the re- 
sponse R shows a maximum as a function of the number 
of units, that follows the same pattern as the maximum M 
(K — 0.2) (inset). Since N decreases the influence of a single 
neighbor, and C increases it, when the coupling intensity is 
larger, the optimal system size increases. Values of T and A 
as in Fig. [3] 



of a larger system is best amplified at a higher probability 
of repulsive links. As the fraction of repulsive links must 
not exceed the fraction of positive ones, there can be a 
limit on how large can a system be, to be able to amplify 
a signal. The same behavior, focusing on the number of 
neighbors was found in a previous study of an Ising-like 
network model [28j . 

IV. SPECTRAL ANALYSIS 

We have already commented that the optimal prob- 
ability of repulsive links drives the system to a glassy 
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FIG. 9: System size induced resonance, as revealed by the 
resonant trajectory at optimal N = 250. Values of T, A and 
C as in Fig. H 



phase. Anderson [§3, HE} has proposed a connection be- 
tween a glass and a delocalization-localization transition, 
relating the existence of many metastable states with a 
localization of modes. From this proposal, we retain the 
idea to work in the eigenspace of the interaction matrix, 
and to look for the fraction of repulsive links where mode 
localization becomes significant. This approach has the 
virtue of not only identifying the steady states, but also 
to shed light onto how the reaction to perturbations is 
sustained and spreads along the system, depending on 
the fraction of repulsive links. In this manner, we hope 
to locate the region where multistability is expected, and 
also to understand the mechanism of response to external 
perturbations. 

Following [27]], let us define the eigenvalues Q a and 
(normalized) eig envectors e a = (e", . . . , e^) of the Lapla- 
cian matrix [29( Jy : 

N 

J'ij = Jij ~ &ij Jkj, (5) 
k=l 

N 

Y / J^ = Q*<- (6) 

3=1 

The effect of the competitive interactions can be de- 
scribed by the so-called participation ratio of eigenvec- 
tor e Q , defined as PR Q = ^ / J2iLi[ e ?] 4 ■ ^ quantifies 
the number of components that participate significantly 
in each eigenvector. A state a with equal components 
has PR Q = N, and one with only one component has 
PR a = 1. When PR a = 1 on a fraction / of elements, 
and elsewhere, then PR a = /, which justifies its name. 
More precisely, we will define "localized" modes as the 
ones whose participation ratio is less than O.liV. Our 
first observation (Fig. [TU)) is that at the optimal region 
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FIG. 10: The participation ratio PR a is a measure of local- 
ization; it estimates the number of eigenvectors components 
contributing to the Q a eigenvalue. The eigenvalues at both 
ends of the spectrum are localized for intermediate levels of 
disorder. N = 200, results after 1000 independent runs. 



p there is a significant fraction of positive eigenvalues, 
and, of those, a significant fraction of the corresponding 
eigenstates are localized. In this region, we will neglect 
the coupling between different modes. This approxima- 
tion allows us to look in more detail at what happens at 
the optimal region, and in particular at the effect of the 
coupling strength C and the number N of elements. 

Let us focus first on the unforced system (A = 0), 
to see how the presence of the disorder induced by the 
repulsive links affects a state configuration. We assume 
each unit i is initially at a given state s°, chosen from a 
random symmetric distribution and split the variables in 
the steady state as Sj = s° + being Xi the deviation 
from the initial condition. We express Xi in the eigenbasis 
of the J/„ matrix: 



N 



a=l 



(7) 



Expanding Eq. (J7XJ) f° r A = 0, multiplying the resulting 
equation by ef, summing over all elements i, and approx- 
imating averages of the product of initial conditions and 
eigenvectors by the product of their individual averages 



(e.g- Eili s ° e 



0), we obtain: 



K - C— 
N 



B a = 0, (8) 



where 



N 



;.=i 



JV 



(9) 
(10) 



1=1 



Neglecting coupling between modes leads to F^ lva = 
1/PR Q if a = p = 7 = 7} and F^ a = otherwise. We 



then obtain the following equation for the amplitude of 
the cv-th mode: 



B. ' + PR a { K - C^j B a = 0. (11) 

According to this approximation, unless Q a > ^ff-, the 
amplitude B a of the mode a is zero, and any small per- 
turbation vanishes. Otherwise, the mode is said "open" 
and B a takes one of the values: 



4, =±,;PR„ [C^ 



K 



(12) 



For intermediate amounts of disorder, some open modes 
begin to appear. The final state of an unit is Si — 
s° + J2a=i B a ef, and when the initial conditions are 
random and the open modes a are localized, the system 
reaches many metastable states, given all the possible 
combinations of individual states. For this reason, we 
want to locate a transition to a region with a significant 
number of localized modes. 

To concretize, we define a measure M of localization: 



M = 



N Q N' 



(13) 



where No is the number of modes a whose associated 
eigenvalue Q a is greater than ^j-, and Nl is the num- 
ber of those modes which, in addition, are localized, i.e. 
PR Q < 0AN. 

Recalling the definition of K (Eq. (TTU10 . we see its value 
is related to a choice of initial conditions, by the variance 
of s°. Since we expect multistability to emerge when the 
initial distribution is more or less uniform, we present 
results in Fig. [TT] for values of K ss 0. 

At moderate levels of disorder, the localized nodes ap- 
pear on the tails of the spectra, Fig. [TO] We confirm that 
the optimal probability of repulsive links coincides with 
a maximal localization of open modes in that region, as 
identified by the peak in M (Fig. fTTjl . 

In a particular metastable state, the units are ran- 
domly distributed, more concentrated near one of the po- 
tential wells. Observing the results in Fig. [TT] for K > 0, 
we see that we recover the dependance on C, and that 
the peak in M still coincides with the optimal probabil- 
ity region. The enhanced responsiveness to an external 
signal can thus be understood as a consequence of mode 
localization. Since units can be in different positions, 
some will be able to answer the signal, and then - since 
the overall coupling is attractive - pull the others. This 
is done in an incremental fashion, as confirmed by the 
localized reaction to perturbations. 

The same analysis is valid when we plot M as a func- 
tion of C or N. We notice a peak in M and accordingly 
the dependance of the response on C (Fig. [6]) and N 
(Fig.[SJ) shows a maximum for intermediate values (insets 
Fig. |5]and Fig. [5]). When C is small, even if the modes 
are open, their amplitude B a is weak, Eq. (fT2"j). A high 
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FIG. 11: Measure of localization M: close to the optimal 
region there is an increase in the fraction of localized open 
eigenstates. The small dependance of M with K and C is as 
expected. Parameters: TV = 200. Inset: A = 0.2, T = 300. 



fraction of repulsive links, increasing the number of open 
modes, can overcome this situation to a certain degree, 
allowing for resonances at a smaller coupling strength. 

V. CONCLUSIONS 

In this work, we have analyzed the response to a weak 
period signal, of a model composed by bistable units cou- 
pled through both attractive and repulsive links. 

Our main result is that the system collective response 
is enhanced by the presence of an intermediate fraction 
of repulsive links. Hence, competitive interactions are 
taken as a source of disorder, as an alternative to pre- 
vious studies where disorder was induced by noise [J] or 
diversity (2lj . and a similar amplification was verified. 

We have chosen a very generic double- well model, and 
have shown that the optimal disorder is the one that 
destroys the ordered system bistability. The resulting 
multistable effective potential allows for the amplifica- 
tion of very weak or fast signals. There is not a need to 
match specific levels of disorder with specific frequencies. 



Having the optimal disorder, the system becomes more 
sensitive to external signals of every kind. Furthermore, 
we have shown that varying the number of elements or 
the coupling strength in an ensemble of coupled bistable 
elements can improve the sensitivity to an external forc- 
ing. These various ways to increase sensitivity make the 
phenomenon less dependent on a fine tuning of the pro- 
portion of repulsive links, which can be a positive feature 
in practical applications. Apparently, when the system 
size becomes very large, it is difficult to get a resonance 
effect, unless we increase the coupling strength by many 
times. Arguably, thi s difficulty can be overcome by other 
types of network settings [361 ] . 

Finally, we have shown that the location of the reso- 
nance peaks can be predicted by a spectral analysis of 
the Laplacian matrix. In heuristic terms [29(, the posi- 
tive eigenvalues of the Laplacian can be seen to express 
the contribution of the coupling term to the vulnerabil- 
ity of the system to perturbations. We conclude that the 
location of the amplification region, for a given system 
size and coupling constant, is reasonably independent of 
the particular dynamical system. In broad terms, it cor- 
responds to the point where the positive eigenvalues of 
the Laplacian matrix become localized, signaling a tran- 
sition to a region where perturbations can accumulate in 
an incremental manner. The more precise location would 
depend on the particular dynamical system by means of 
a condition on open modes. 

Competitive interactions are widespread in nature, no- 
tably in biological systems. In those systems and oth- 
ers, there has been some studies highlighting their role 
in achieving a coherent behavior in the absence of forc- 
ing: increa sing synchronization (37l | or enabling a collec- 
tive firing [30| . In the present study, we saw they can 
also help to enhance perception, something that can be 
potentially relevant in sensory systems. 
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